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Summary. The transport of products between different suppliers or production 
units can be described similarly to driven many-particle and traffic systems. We 
introduce equations for the flow of goods in supply networks and the adaptation of 
production speeds. Moreover, we present two examples: The case of linear (sequen- 
tial) supply chains and the case of re-entrant production. In particular, we discuss 
the stability conditions, dynamic solutions, and resonance phenomena causing the 
frequently observed "bullwhip effect", which is an analogue of stop-and-go traffic. 
Finally, we show how to treat discrete units and cycle times, which can be applied 
to the description of vehicle queues and travel times in freeway networks. 
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Supply chain management is a major subject in economics, as it signif- 
icantly determines the efficiency of production processes [1]. Many related 
studies focus on subjects like optimum buffer sizes and stock levels, but a 
stable dynamics and optimal network structure are also important subjects 
[2,3]. Therefore, this scientific field is connected with the statistical physics 
of networks [4]. 

In this paper, we will investigate the dynamic properties and linear stabil- 
ity of supply chains by a generalization of "fluid-dynamic" production models, 
which have been inspired by traffic models [5-10]. Fluid-dynamic models take 
into account non-linear interactions and are suitable for on-line control, as 
they are numerically much more efficient than event-driven (Monte-Carlo) 
simulations. Moreover, unlike queueing theory, they are not mainly restricted 
to the treatment of stationary situations, but can reflect variations in the con- 
sumption rate and effects of machine breakdowns or changes in the production 
schedule. 
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The organization of this contribution is as follows: Section 1 focusses on 
the description of supply and production systems, in particular linear supply 
chains and re-entrant production. Section 2 treats freeway traffic as a queueing 
network and applies the formula for the cycle time to the determination of 
travel times. 



1 Modeling production systems as supply networks 

Our production model assumes u production units or suppliers j which deliver 
dij products of kind i G {1, . . . ,p} per production cycle to other suppliers and 
consume Ckj goods of kind k per production cycle. The coefficients Ckj and 
dij are determined by the respective production process, and the number of 
production cycles per unit time (e.g. per day) is given by the production 
speed Qj(t). That is, supplier j requires an average time interval of 1/Qj(t) 
to produce or deliver d^ units of good i. The temporal change in the number 
Ni(t) of goods of kind i available in the system is given by the difference 
between the inflow 



Qt(t) = J2 d Mt) (1) 



and the outflow 



u+l 



g° ut (*) = E c ^(*)- ( 2 ) 

In other words, it is determined by the overall production rates dijQj(t) of 
all suppliers j minus their overall consumption rates CijQj(t): 

= Qt(t) - Q° ut (t) = X>, - dj)Qj{t) Yi(t) . (3) 

Herein, the quantity 

Yi(t) = Cuu+iQu+i{t\ - d iQ Qo{t) (4) 

consumption and losses inflow of resources 

comprises the consumption rate of goods i, losses, and waste (the "export" of 
material), minus the inflows into the considered system (the "imports"). In 
the following, we will assume that the quantities are measured in a way that 
< Cij , dij < 1 (for 1 < i < p, 1 < j < u) and the "normalization conditions" 

u u 

djo = l- y~] d i:j > , c iiU+ i = 1 - ~y] Cjj > (5) 
i=i j=i 

are fulfilled. Equations (3) can then be interpreted as conservation equations 
for the flows of goods. 
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1.1 Adaptation of production speeds 

The production speeds Qj(t) may be changed in time in order to adapt to 
varying consumption rates Yi(t). We will assume that it takes a typical time 
period T to adapt the actual production speeds Qj(t) to the desired produc- 
tion ones, Wj\ 



Here, the curly brackets indicate that the so-called management or control 
function Wj(. . . ) may depend on all inventories Ni(t) with i G {1, . . . ,p}, their 
derivatives dNi/dt, and/or all production speeds Qk(t) with k G {1, . . . ,u}. 
The resulting dynamics of the supply network can be investigated by means of 
a linear stability analysis, which has been carried out for linear supply chains 
in Ref . [2] and for supply networks with dij = 5ij in Rcf . [3] . Some of the main 
results will be discussed in the following. 

1.2 Modelling one-dimensional supply chains 

For simplicity, let us investigate a model of onc-dimensional supply chains, 
here, which corresponds to dij = % and Cjj = where Sij = 1 for i = j 

and % = otherwise. This implies Qf(t) = Qi (t] I and Q° ut (t) = Q»+i(*) = 
Qf\_i(t). The assumed model consists of a series of u suppliers i, which receive 
products of kind i — 1 from the next "upstream" supplier i — 1 and generate 
products of kind i for the next "downstream" supplier i + 1 at a rate Qi (t) 
[7,9]. The final products are delivered at the rate Q u {t) and removed from 
the system with the consumption rate Y u (t) = Q u+ i(t). The consumption 
rate is typically subject to perturbations, which may cause a "bullwhip effect" 
[11-18], i.e. growing variations in the stock levels and deliveries of upstream 
suppliers. This is due to delays in the adaptation of their delivery rates. 
The inventory of goods of kind i changes in time t according to 



^ = ^[w j ({N i (t)},{dN i /dt},{Q k (t)})-Q j (t) . 



(6) 



_i = Qf{t) Q™\t) = Qi{t) Q i+1 {t) , 



(7) 



while the temporal change of the delivery rate will be assumed as 



^ = ^[Wi(Ni(t), dNi/dt, Qi(t),Q i+1 



(t))-Qi(t) 



(8) 



The management strategy 




^+e[Q°i-Qm) 



(9) 



appears to be appropriate to keep the inventories Ni (t) stationary, to maintain 
a certain optimal inventory Nf (in order to cope with stochastic variations due 
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to machine breakdowns etc.), and to operate with the equilibrium production 
rates Q° = Y®, where denotes the average consumption rate. The supplier- 
independent parameter values r, (3, and e can be justified with suitable scaling 
arguments [3]. 



1.3 Dynamic solution and resonance effects 

In the vicinity of the stationary state characterized by TV, (t) — N® and Qi (t) — 
it is possible to calculate the dynamic solution of the one-dimensional 
supply chain model [3, 8]. For this, let rii{t) — Ni{t) — Nf be the deviation of 
the inventory from the stationary one, and q%(t) = Qi(t) — the deviation 
of the delivery rate. The linearized model equations read 



and 



dq% 1 / rii drii 



(10) 



(11) 



Deriving Eq. (11) with respect to t and inserting Eq. (10) results in the fol- 
lowing set of second-order differential equations: 



d 2 qi (0+j) dfH J_ , , _ 1_ 
dt 2 + T dt + Tt qA ) ~ T 



dt 



(12) 



=2 7 



This corresponds to the differential equation for the damped harmonic oscil- 
lator with damping constant 7, eigenfrequency u, and driving term fi(t). The 
two eigenvalues of this system of equations are 



(13) 



i.e. for (fi + e) > their real parts are always negative, corresponding to a 
stable behavior in time. Nevertheless, we will identify a convective instability 
below, i.e. the oscillation amplitude can grow from one supplier to the next 
one upstream. 

The set of equations (12) can be solved successively, starting with i = u and 
progressing to lower values of i. For example, assuming periodic oscillations 
of the form f u {t) = /° cos(crf), after a transient time much longer than I/7 
we find 

q u (t) = f°Fcos(at + <p) (14) 



with 



tamp 



27 a 



a(P + e) 
a 2 T- 1/t 



(15) 
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and 

1 T 

F = , = , , (16) 

^{a 2 - lu 2 ) 2 + 4 7 2 a 2 \J\a 2 T — 1/r] 2 + a 2 {(3 + e) 2 

where the dependence on the eigenfrcquency u> is important to understand 
the occuring resonance effect. Equations (12) and (14) imply 



fu-l(t) — — 



quit) R dq u 
r +l dt 



/°_ lC0S (ai + ^ + <5) (17) 



with 



tanS = af3r and /„% = /°- V(l/r) 2 + (a/3) 2 . (18) 



1.4 "Bull- whip effect" and stop-and-go traffic 



The oscillation amplitude increases from one supplier to the next upstream 
one, if 

_ L , a 2 K £ + 2/3)-2r/r] + ct 4 r 2 |- 1 / 2 

/2 ~ I (Vr) 2 + M) 2 J >1 - (W) 

One can see that this resonance effect can occur for < a 2 < 2/(Tr) — e(e + 
2f})/T 2 . Therefore, variations in the consumption rate are magnified under 
the instability condition 

T > €T (J3 + e/2) . (20) 

Supply chains show this "bullwhip effect" (which corresponds to the phe- 
nomenon of convective, i.e. upstream moving instability), if the adaptation 
time T is too large, if there is no adaptation to the equilibrium production 
speed j corresponding to e — 0, or if the management reacts too strong to 
deviations of the actual stock level Ni from the desired one N®, correspond- 
ing to a small value of r. The latter is very surprising, as it implies that the 
strategy 

= ±[N? - NM, (21) 

which tries to maintain a constant work in progress Ni(t) — N®, would ulti- 
mately lead to an undesireable bullwhip effect. In contrast, the management 
strategy 

f -i{-^ + «Wf -*(.)!} <»> 

would avoid this problem, but it would not maintain a constant work in 
progress. The control strategy (9) with a sufficiently large value of r would 
fulfill both requirements. 

The "bullwhip effect" has, for example, been reported for beer distribution 
[22, 23], but similar dynamical effects are also known for other distribution or 
transportation chains with significant adaptation times T. It has, for example, 
some analogy to stop-and-go traffic [3, 10], where delayed adaptation also leads 
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to an unstable behavior. In the case (3 — and e = 1, the stability condition 
(20) agrees exactly with the one of the optimal velocity model [24], which is 
a particular microscopic traffic model. This car-following model assumes an 
acceleration equation of the form 

dvi(t) _ V op M(t)) - Vl (t) (23) 



dt 

and the complementary equation 

ddi(t) 



dt 



= -[vi{t)-v i+1 {t)}. (24) 



In contrast to the above supply chain model, however, the index i represents 
single vehicles, Vi(t) is their actual velocity of motion, V opt the so-called op- 
timal (safe) velocity, which depends on the distance di(t) to the next vehicle 
ahead. T denotes again an adaptation time. Comparing this equation with 
Eq. (9), the velocities Vi would correspond to the delivery rates Qi, the opti- 
mal velocity V op t to the desired delivery rate W, (Ni) = Q® + (Nf — Ni) / t, and 
the inverse vehicle distance 1/di would approximately correspond to the stock 
level Ni (apart from a proportionality factor). This shows that the analogy 
between supply chain and traffic models concerns only their mathematical 
structure, but not their interpretation, although both relate to transport pro- 
cesses. Nevertheless, this mathematical relationship can give us hints, how 
methods, which have been successfully applied to the investigation of traffic 
models before, can be generalized for the study of supply networks. Compared 
to traffic dynamics, supply networks and production systems have some inter- 
esting new features: Instead of a continuous space, we have discrete production 
units j, and the management strategy (6) is generally different from the ve- 
locity adaptation (23) in traffic. With suitable strategies, in particular with 
large values of r and /3, the oscillations can be mitigated or even suppressed. 
Moreover, production systems are frequently supply networks with complex 
topologies rather than one-dimensional supply chains, i.e. they have additional 
features compared to (more or less) one-dimensional freeway traffic. They are 
more comparable to street networks of cities [25]. 



1.5 Calculation of the cycle times 

Apart from the productivity or throughput Qi of a production unit, produc- 
tion managers are highly interested in the cycle time Ti, i.e. the time interval 
between the beginning of the generation of a product and its completion. Let 
us assume that the queue length Li (t) of products waiting to be processed by 
production unit i is given by the inventory Ni-i(t) of product The change 
of the queue length Li (t) in time is then determined by the difference between 
the arrival rate Qf TT (t), which corresponds to the inflow Q™_i = Qi-i(t) from 
production unit i — 1, and the departure rate Qi° p (t), which corresponds to 
the outflow Q°"\(t) — Qi(t) to production unit i: 
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dJ i = dJ ir = Qi - l{t) Qi{t) = g ™ Q « CP(<) ■ (25) 

On the other hand, the waiting products move forward Qf ep (t) = Qi(t) steps 
per unit time, as Qi(t) is the processing rate (production speed). For this 
reason, the overall time Ti(t) until having been processed is given by the 
implicit equation 

t+T,(t) t+T,(t) t 

JVi-i(t) - Ht) = J dt' of p (0 = J dt' of p (t') - | dt' of p (0 , 

i — oo — oo 

(26) 

if the queue of length Li{t) was joined at time t. Accordingly, high inventories 
imply long cycle times, which favours just-in-time production with small stock 
levels. From Eqs. (25) and (26), one can finally derive a delay-differential 
equation for the waiting time under varying production conditions [10]: 

dJ\ = Qf"(t) _ = Q a -!(t) _ 
dt Q^(t + Ti(t)) Qi(t + Ti(t)) ■ 1 ) 

This equation can be solved numerically as a function of the production 
rates Qi(t), since the production initially starts with a cycle time of Tj(0) = 
l/0i(0), corresponding to the average processing time 1/Qj(0) when the pro- 
duction unit i is started. In this way, it is possible to determine the process 
cycle times T and the overall production time as the sum of the cycle times 
of all single production steps, taking into account the respective time delays 
Tf. When a specific product enters the queue before production unit i at time 
t = ti, it enters the queue before production unit i+1 at time ti+i = ti+Ti(U). 
The time of delivery to the customer is given by t u +i = t u + T u (t u ). 



1.6 Modeling of discrete units 

The above "fluid-dynamic" model equations can not only be used to repre- 
sent approximate mean values of large numbers of products. They can also be 
transfered to the treatment of discrete units (such as single units of a prod- 
uct), if their dynamics is sufficiently deterministic. For example, one could 
represent the time interval AT, during which a discrete unit occupies a cer- 
tain production unit i by a step function. However, as step functions are not 
diffcrentiable everywhere, we will replace them by smooth functions. 

One possible specification uses a Fourier approximation of the step func- 
tions. To lowest order, one may take 

Qi(t) - A]T{1 - cos^t - t k )/AT]} (28) 
k 

where t\ denotes the starting time of occupation by object k. The cosinus 
function is set to zero for t <t\ and t > t\ + AT. The prefactor A is deter- 
mined in a way that satisfies the normalization condition, 
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t'i+AT 

J dt'Qi(t') = l, (29) 

which implies A — 1/AT. Under this condition, the integral 

N™(t) = / dt'Qdt 1 ) (30) 
Jo 

counts the number of units that have passed the cross section or production 
unit under consideration. 

Another possible specification assumes 

Q t (t)=BY / (t-t>;) 2 (t-t*-AT)\ (31) 

k 

where we set the fourth order polynomial to zero for t < t\ and t > t\ + AT. 
The normalization condition (29) implies B = 30/(Z\T) 5 . 

Together with (28) or (31), the relationship (30) can also be applied to 
situations where several units occupy a production unit at the same time. As 
before, the management function can be chosen as a function of the inventories 
Ni(t), but in many cases, it would be reasonable to replace a reaction to 
temporal changes dNi/dt in the inventories by exponentially smoothed values. 



1.7 Re-entrant production 



Clin 



Qout r~\a.r 
VI 



Re-entrant 
Production 

N 



Q'2 



Ni 



QT 1 



N 2 



-\ov,t 
{2 



Qout 
u-1 



QT P Ql n 



Fig. 1. Illustration of the model of re-entrant production and its variables (see main 
text for details). 

Semiconductor production [19] relies on some extremely expensive (litho- 
graphic) production units, which are therefore used for many similar subse- 
quent production steps i, namely the production of the various layers of a 
chip. This is the reason for re-entrant production [20,21] (see Fig. 1.7), pos- 
ing particular control problems not only for the overall arrival rate Q arr (£) 
in the re-entrant production area, but also concerning the fractions Pi(t) of 
products in the different production stages i fed into it. If iV»(i) denotes the 
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stock level of chips after the ith entry (production step), Q™{t) the respective 
inflow, and Q° nt (t) the outflow, the related balance equation is again 

dN- 

-^ = QT(t)-Qr\t). (32) 

Assuming also a buffer before the first and a buffer after the last (wth) produc- 
tion step, this equation applies to i £ {0,1, ... ,u}, where Q° ut corresponds to 
the removal rate of products that have completed the re-entrant (lithographic) 
production steps. The outflow Q°H\(t) determines the arrival rate Qf"(t) for 
the ith re-entrant production step, while the departure rate Q^ cp (t) deter- 
mines the inflow Q™{t) into the subsequent buffer. The overall stock level 
N(t) of chips at various production stages in the re-entrant production area 
changes in time according to 

rlN U U 

= £[«(*) - = £[<?™ - Qt cp (t)\ = Q aiI (t) Q dcp (t) ■ 

i=l i=l 

(33) 

The overall departure rate Q dcp (t) from the re-entrant production area de- 
pends on the number N(t) of products processed at the same time. Moreover, 
if pi (t) denotes the fraction of chips that have entered the re-entrant area for 
the ith production step at time t and T° (t) is the related cycle time of these 
chips in the re-entrant area, their departure rate is given by 

Qf p (t) = Pi (t - T"(t))Q d ^(N(t)) = Qf(t) . (34) 

As before, the delay-differential equation for the temporal change of the cycle 
time is 

dt Qtep(t + T°(t)) ' 1 ' 

The overall arrival rate may be adapted according to 
r/O arr 1 

= j;[W(N(t),dN/dt,Q^(t),Q d ^(t),N u (t),dN u /dt,Q^(t))-Q a "(t)} , 

(36) 

where the management function W(. . . ) depends not only on the variables 
N(t) and Q dcp (t) characterizing the re-entrant production area, but also on 
the (desired) removal rate Q° ut (<) and the stock level N u (t) of the final buffer. 
Finally, the specific arrival rates of chips for the (i+l)st re-entrant production 
step are 

Qf^(t) = Pi +1 (t)Q m (t) = Q° Ut (t) , (37) 
where their relative percentages 

(f] = Wi (Nj (t) , dNj/dt, Q?(t),Q™\t), Q° ut (t) ) 
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are controlled by another management function Wj(. . . ), which depends on 
the stock levels iVj(i) of the respectively preceding buffers and on the removal 
rate Q° ut (t). The management functions Wi{. . .) allow one to specify differ- 
ent priorities such as push or pull strategies. Specifications of the management 
function, which can avoid bullwhip (resonance) effects and long cycle times 
will be presented in a forth-coming paper. Generalizations to the simultane- 
ous production of various products are possible, requiring the consideration of 
an individual number of re-entrant steps with potentially different production 
speeds and cycle times ("overtaking") in the re-entrant area. The problem is 
similar to the treatment of heterogeneous multi-lane traffic [26] . In the follow- 
ing, we will present a simpler traffic model for uniform vehicles, which shows 
how to transfer the above approach from production processes to transporta- 
tion processes in street networks. 



2 Queueing model of vehicle traffic in freeway networks 

In the following, we will propose a traffic model [25], which was inspired by 
the above model of supply networks. Here, the elements to be served are the 
vehicles. The formula for the queue length determines the density in a road 
section, and the cycle time corresponds to the travel time. 

When we now specify road traffic as a queueing system, we will take into 
account essential traffic characteristics such as the flow-density relation or 
the properties of extended congestion patterns at bottlenecks. In fact, traffic 
congestion is usually triggered by spatial inhomogeneities of the road network 
[27], and queueing effects are normally not observed along sections of low 
capacity, but upstream of the beginning of a bottleneck. Therefore, we will 
subdivide roads into sections i of homogeneous capacity and length Z™ ax , 
which start at place Xi and end with some kind of inhomogeneity (i.e. an 
increase or decrease of capacity) at place Xi + i — Xi + Z*™ ax . In other words, 
the end of a road section i is, for example, determined by the location of an 
on- or off-ramp, a change in the number of lanes, or the beginning or end of 
a gradient. 

The model can be derived from the "fluid-dynamic" continuity equation 

' — + ^MEl ^ — source terms (39) 
at ox 

describing the conservation of the number of vehicles. Here, p(x,t) denotes 
the vehicle density per lane at place x and time t, and Qi(x,t) the traffic 
flow per lane. The source terms originate from ramp flows Ql anip {t), which 
enter the road at place Xi+i. Let us define the arrival rate at the upstream 
end of road section i (the "inflow") by Qf"(t) — hQi{xi + dx,t), where dx 
is a differential space interval and J, the number of lanes of road section i. 
Analogously, the departure rate from the downstream end of this section is 
defined by Q^ cp (t) = IiQi(xi + i~dx, t). (Note that this "outflow" from section 
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i is to be distinguished from the outflow Qout per lane from congested traffic) . 
The conservation of the number of vehicles implies that the departure rate 
plus the ramp flow determine the arrival rate in the next downstream section 
i + l: 

In order to guarantee non-negative flows, we will demand for the ramp flows 
that the consistency condition -Q*° p (t) < Qj amp (t) < Qf^{t) is always met. 

Integrating the continuity equation over x with Xi < x < x i+ i provides a 
conservation equation for the number Ni(t) — J^ +1 dx I%p{x, t) of vehicles in 
road section i. It changes according to 

= Qr{t ) Q ^(t) = Qtm + QtT(t) - Qt p (t) ■ (4i) 

In the terminology of queueing theory, equation (41) reflects the change in the 
number Ni(t) of vehicles waiting to be served by the downstream end of the 
section, and the number Ii of lanes corresponds to the number of "channels" 
serving in parallel. Moreover, iVj(i) corresponds to the "queue length". How- 
ever, this does not necessarily mean that we have congested traffic as, in the 
terminology of traffic theory, "queue length" refers to something else, namely 
the spatial extension ?j(£) > of a congested road section. In the following, 
we will try to express the traffic dynamics and the travel times only through 
the flows at the cross sections x i} taking into account the features of traffic 
flow in a simplified way. 

For free flow, i.e. below some critical vehicle density p cv per lane, the re- 
lation between the traffic flow Qi per lane and the density p per lane can be 
approximated by an increasing linear relationship, while above it, a falling lin- 
ear relationship is consistent with congested flow-density data (in particular, 
if the average time gap T is treated as a time-dependent, fluctuating variable) 
[28]. This implies Qi(x,t) « Q l {p{x,t)) with 



- i Q con g(p) = (1 _ p/p . am)/T otherwise. 



Here, V® denotes the average free velocity, T the average time gap, and pj am 
the density per lane inside of traffic jams. Moreover, we define the free and 
congested densities by 

P? ee (Qi) = Qi/V° and P T K (Qi) = {l-TQ i )p iBja . (43) 

The quantity Q ou t = (1— Pa/ Pjam)/T corresponds to the outflow per lane from 
congested traffic [29]. Depending on the parameter specification, the model 
describes a continuous flow-density relation (for p C rV® = Qout) or a capacity 
drop at the critical density p cr and high-flow states immediately before (if 
PcrV° > Q out ). 

According to shock wave theory [30] , density variations at place x propa- 
gate with velocity C(t) = [Qi{x + dx, t) — Qi(x — dx, t)]/[p(x + dx, t) — p(x — 
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dx, t)]. Accordingly, the propagation velocity is C — V- 3 in free traffic, and 
C = — c = —l/(Tp max ) in congested traffic. Therefore, it takes the time pe- 
riod T/ ree = ^™ ax /^i° f° r a perturbation to travel through free traffic, while 
it takes the time period T- cong = /™ ax /c, when the entire road section i is 
congested. 

Now, remember that congestion in section i starts to form upstream of a 
bottleneck, i.e. at place Xi+i. Let k(t) denote the length of the congested area 
and x(t) — Xi+i — U(t) = Xi + Z' nax — U(t) the location of its upstream front. 
Then, we have free traffic between Xi and Xi + l™ ax — k(t), i.e. Qi(x — dx, t) = 
Qi"(t — ( x — X i)/Vi) (considering dx — > 0), and congested traffic downstream 
ofx(t), i.e. Q^x+dx.t) = Qf p (t-(x l+1 -x)/c). Withdx/dt = -dk/dt = C(t) 
and Eq. (43) we find 

dt P r s (Q- cp (^ - h{t)/ c )/ii) - p?°HQrtt - K max - imm/h) ' 

(44) 

The capacity of a congested road section i is approximated as the outflow 
<3out = (1 — Pci/Pjam)/T from congested traffic per lane times the number li 
of lanes, minus the maximum bottleneck strength at the end of this section. 
This may be given by an on-ramp flow Q^ amp (t) > or analogously by (Ii — 
Ii+i)Qout in case of a reduction I i+ i — Ii < in the number of lanes, or in 
general by some time-dependent value AQi(t) in case of another bottleneck 
such as a gradient: 

QT P (t) - liQout - max[Q: amp (t), {h - I i+1 )Q out ,AQi(t), 0] . (45) 

Analogously, the maximum capacity Qf 1 ax (i) of the road section i under free 
flow conditions is given by the maximum flow Iip CI V° minus the reduction by 
bottleneck effects: 

Qr x (t) - I lPcl V° - max[Qr mP W, {h - h + i)p^, AQ^O] . (46) 

Moreover, we have to specify the departure rate Qi° p (t) as a function 
of the respective traffic situation. Focussing on the cross section at location 
Xi + i and considering the directions of information flow (i.e. the propagation 
direction of density variations), we can distinguish three different cases: 

1. If we have free traffic in the upstream section i and free or partially con- 
gested traffic in the downstream section i + 1, density variations propagate 
downstream and the departure rate Q^ ep (t) at time t is given as the arrival 
rate Qf rT (t - 7f cc ) = Q^ p (t - 7f cc ) + Q™" p (i - T/ ro °), since the vehicles 
entering section i at time t — T- ree leave the section after an average travel 
time T t of lf ee . 

2. In the case of partially or completely congested traffic upstream and free 
or partially congested traffic downstream, the departure rate Q^ cp (t) is 
given by the capacity Qi ap {t) of the congested road section i. 
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3. In the case of congested traffic on the entire downstream road section 
z + l, the departure rate Q^° p (t) is given by the departure rate Qi^{t — 
j,con g ^ £ rom £jj e downstream section at time t — T™" s minus the ramp flow 
Q T i amp (t) entering at location x i+ i. 

Summarizing this, we have 

(Qr(t - T/ rcc ) if l l+1 (t) < and k(t) = 

QT P (t) = { Q? P (t) if h+i{t) < IT+r and k{t) > (47) 

A numerical solution of the above defined section-based qucucing-theore- 
tical traffic model is carried out as follows: First, calculate the new arrival 
and departure rates by means of Eqs. (40) and (47), taking into account 
the boundary conditions for the flows at the open ends of the road network. 
Second, determine the queue lengths k(t) in all road sections i: (i) If traffic in 
the road section flows freely (h(t) = 0) and the maximum capacity <2™ ax (i) 
is not reached, i.e. Qf"(t - Tf rco ) < Q™ ax (t), we have dk(t)/dt = and 
the traffic flow in the road section remains free, (ii) If the road section is 
completely congested (h(t) = i™ ax ) and the arrival rate Qf"(t) is not below 
the departure rate at time t - T; cong , i.e. Qf p (t - T™ ng ) < Qf"(t), the road 
section i stays fully congested and dk/dt = 0. (iii) In other cases, we have 
partially congested traffic in road section i and the length k (t) of the congested 
area changes according to Eq. (44). Next, one continues with the first step 
for the new time t + dt, and so on. It is obvious, that this numerical solution 
is significantly more simple and robust than the numerical solution of the 
Lighthill-Whitham model, as shock waves (i.e. the interfaces between free and 
congested traffic) are treated analytically and the propagation velocities of 
perturbations within the free and congested regions are constant. 

The travel time Ti(i) of a vehicle that enters road section i at time t can be 
calculated analogously to Eq. (27) [25], i.e. via the delay-differential equation 

dm) = Qr(t) 1 = otm+Q^jt) i 

dt " Q^(t + Ti(t)) " Qt P (t + m)) 

According to this, the travel time Ti(t) increases with time, when the arrival 
rate Qf" at the time t of entry exceeds the departure rate Q^ ep at the leaving 
time t + Ti(t), while it decreases when it is lower. It is remarkable that this 
formula does not explicitly depend on the velocities on the road section, but 
only on the arrival and departure rates. The calculation of the travel time 
based on the velocity is considerably more complicated: Let v(t) be the velocity 
and x(t) = Xi + J* dt' v(t') the location of a vehicle at time t, when it 
enters section i at time to- Its travel time Tj(io) on section i is given by the 
implicit equation x(t + Ti(t )) = x i+ll which says that this vehicle reaches 
place Xi+i = Xi + Z*™ ax at time to + Ti(to). The vehicle speed v(t) is also 
difficult to determine, as it depends on the (free or congested) traffic state at 
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its respective location x(t): It is v(t) — V® in free flow, i.e. for V°(t — t ) < 
Z 4 max - k(t). In congested flow, i.e. for t + ^ max ~ k(t)]/V? < t < i +T 4 (i ), it 
is determined viav(t) = Qi(x(t),t)/ p(x(t),t) — [l/p(x(t), t) — \/p max \/T with 
p(x(t),t) = p(x i+1 ,t-[x i+1 -x(t)]/c) = {l-TQf cp (t-[x i+1 -x(t)]/c) I Ii}p i!im . 

3 Summary and conclusions 

In this contribution, traffic and production systems have been treated in a 
uniform way as dynamic queueing networks, since it does not matter whether 
one treats the transport of goods from one production unit to the next one 
or of vehicles from one cross section of the road network to the next. Conse- 
quently, the stability conditions for supply chains looked similar to the ones of 
some specific traffic models. Moreover, we have presented a delay-differential 
equation for the determination of cycle (production) times, which can be also 
applied to calculate the travel times of vehicles. Interestingly enough, this 
formula did not require to calculate the vehicle speed on freeway sections, 
but only the in- and outflows at certain cross sections of the street network, 
namely, where the road capacity changed. 

The derived instability conditions allow one to choose appropriate manage- 
ment strategies which can avoid the well-known bullwhip effect. This effect 
describes an amplification of variations in the delivery rate and inventory 
from one supplier to the next one upstream. Such a convective instability can 
occur despite of a stable behavior in time because of the possibility of reso- 
nance effects. Apart from this, we have sketched the treatment of re-entrant 
production processes and of discrete units. 

The advantage of "fluid-dynamic" models of traffic, supply, and produc- 
tion networks is their great numerical efficiency and their consideration of 
non-linear interaction effects. Therefore, in contrast to most classical queue- 
ing theoretical approaches and to event-driven (Monte Carlo) simulations, 
they are suitable for on-line control and the treatment of variations in the 
consumption rate or in the production program. 

Present research focusses on the effect of the topology of supply networks 
on their dynamics [3]. The dynamics of production processes can even be 
chaotic [22,31-33]. Other studies concentrate on the subject of optimal con- 
trol (including chaos control) [33-39] , which is particularly challenging for re- 
entrant production [20,21]. Finally, the presented traffic model is now being 
applied to the simulation of city networks with adaptive traffic light control 
and to dynamic assignment problems. 
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